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Abstract 

In [5], Srijimtongsiri and Vavasis propose the Kantorovich-Test 
Subdivision algorithm, or KTS, which is an algorithm for finding all 
zeros of a polynomial system in a bounded region of the plane. This 
algorithm can be used to find the intersections between a line and a 
surface. The main features of KTS are that it can operate on poly- 
nomials represented in any basis that satisfies certain conditions and 
that its efficiency has an upper bound that depends only on the con- 
ditioning of the problem and the choice of the basis representing the 
polynomial system. 

This article explores in detail the dependence of the efficiency of 
the KTS algorithm on the choice of basis. Three bases are considered: 
the power, the Bernstein, and the Chebyshev bases. These three bases 
satisfy the basis properties required by KTS. Theoretically, Chebyshev 
case has the smallest upper bound on its running time. The compu- 
tational results, however, do not show that Chebyshev case performs 
better than the other two. 
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1 The line-surface intersection problem and 
the required basis properties 

Let (f>o, . . . , <f) n denote a basis for the set of univariate polynomials of degree 
at most n. For example, the power basis is defined by 4>i(t) = t l . The line- 
surface intersection problem can be reduced to the problem of finding all 
zeros of 

m n 

f(u,v) = ^^Cijfaiu)^), 0<U,V<1, (1) 
i=0 j=0 

where Cij GK 2 (i — 0, 1, . . . , m;j — 0, 1, . . . , n) denote the coefficients [5]. 

For this article, let the notation ||-|| refer specifically to infinity norm. 
Other norms are explicitly notated so. The Kantorovich-Test Subdivision 
algorithm (KTS in short), proposed by Srijuntongsiri and Vavasis [5], can be 
used to solve ([I]). KTS works with any polynomial basis (f>i(u)4>j(v) provided 
that the following properties hold: 

1. There is a natural interval [I, h] that is the domain for the polynomial. 
In the case of Bernstein polynomials, this is [0, 1], and in the case of 
power and Chebyshev polynomials, this is [—1, 1]. 

2. It is possible to compute a bounding polytope P of S = {f(u,v) : / < 

u,v < h}, where f(u,v) = YhLo X^jLo c ij ( l ) i{ u ) ( t ) j( v ) and c n e ^ d for 
any d > 1, that satisfies the following properties: 

(a) Determining whether G P can be done efficiently (ideally in 
0(mn) operations). 

(b) The polytope P is affinely invariant. In other words, the bounding 
polytope of {Af(u, v) + b : I < u, v < h} is {Ax + b : x 6 P} for 
any nonsingular matrix A G M. dxd and any vector b G M. d . 

(c) For any y G P, 

\\y\\<6 max \\f(u,v)\\, (2) 

l<u,v<h 

where 9 is a function of m and n. 

(d) If d = 1, then the endpoints of P can be computed efficiently 
(ideally in 0(mn) time). 
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3. It is possible to reparametrize with [I, h} 2 the surface Si = {f(x) : x G 
B(x°,r)}, where x° G M? and r G K > 0. In other words, it is possible 
(and efficient) to compute the polynomial / represented in the same 
basis such that S\ = {f(x) : x G [I, h] 2 }. 

4. Constant polynomials are easy to represent. 

5. Derivatives of polynomials are easy to determine in the same basis, 
(preferably in 0(mn) operations). 

We are generally interested in the case where d = 2. In this case, we call P 
a bounding polygon. Recall that P is a bounding polygon of S if and only if 
x G S implies x G P. 

2 The Kantorovich-Test Subdivision algo- 
rithm 

The description of KTS, as well as the definitions of the quantities mentioned 
in the description, are given below. More details can be found in [5]. 

For a given zero x* of polynomial /, let uj*(x*) and p*(x*) be quantities 
satisfying the conditions that, first, uj*(x*) is the smallest Lipschitz constant 
for fix*)- 1 /', i-e., 

H/Vr 1 (/'(*) " f'(v))\\ < ■ \\x - V\\ for all x, y G 5{x*,p.(x*)) 

(3) 

and, second, 
Define 



7 (0) = 1/ (4^/0(40 + 1) - 80 



where is as in (j2j). Define ujd 1 to be the smallest nonnegative constant w 
satisfying 

\\f\x*)-\f\y)-f\z))\\ < u-\\y-z\\, y, z G D' , x* G [0, l] 2 

satisfying /(x*) = 0, 

where 

/y=[- 7 (0),l + 7(0)] 2 . (6) 
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Denote uuf as the maximum of ud> and all uj^x*) 



ujf = max{co>n/, max ujJx*)}. 

1 x*€C 2 :f(x*)=0 

Finally, define the condition number of / to be 

cond( f) = max{w f , max II f'fx*) -1 f'(y) II}. (7) 

We define the Kantorovich test on a region X = B(x°,r) as the appli- 
cation of Kantorovich's Theorem on the point x° using B(x ,2-f(6)r) as the 
domain (refer to [H H] for the statement of Kantorovich's Theorem). The 
region X passes the Kantorovich test if t]uj < 1/4 and B(x°,p-) C D'. 

The other test KTS uses is the exclusion test. For a given region X, let 
fx be the polynomial in the basis 4>i(u)(pj(v) that reparametrizes with [/, h] 2 
the surface defined by / over X. The region X passes the exclusion test if 
the bounding polygon of {fx(u,v) : I <u,v < h} excludes the origin. 

Having defined the above prerequisites, the description of KTS can now 
be given. 

Algorithm KTS: 

• Let Q be a queue with [0, l] 2 as its only entry. Set S = 0. 

• Repeat until Q = 

1. Let X be the patch at the front of Q. Remove X from Q. 

2. If X £ X s for all X s e 5, 

— Perform the exclusion test on X = B(x°,r) 

— If X fails the exclusion test, 

(a) Perform the Kantorovich test on X 

(b) If X passes the Kantorovich test, 

i. Perform Newton's method starting from x° to find a 
zero x* . 

ii. If x* Xs for any Xs G S (i.e., x* has not been found 
previously) , 

* Compute p*(x*) and its associated u*(x*) by binary 
search. 
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* Set S = Su{B(x*,p*(x*))}. 

(c) Subdivide X along both u and t>-axes into four equal sub- 
regions. Add these subregions to the end of Q. 



The following theorem shows that the efficiency of KTS has an upper 
bound that depends only on the conditioning of the problem and the choice 
of the basis. 

Theorem 2.1. Let f(x) = f(u, v) be a polynomial system in basis <fti(u)<f)j(v) 
in two dimensions with generic coefficients whose zeros are sought. Let 
X = B(x°,r) be a patch under consideration during the course of the KTS 
algorithm. The algorithm does not need to subdivide X if 



Proof. See 0. □ 

Remark: Both terms in the bound on the right-hand side of ([S]) are 
increasing as a function of 1/6. Therefore, our theorem predicts that the 
KTS algorithm will be more efficient for 6 as small as possible (close to 1). 

3 Properties of the power, Bernstein, and 
Chebyshev bases 

As mentioned above, the basis used to represent the polynomial system must 
satisfy the properties listed in Section [I] for KTS to work efficiently. Three 
bases, the power, Bernstein, and Chebyshev bases are examined in detail. 
The power basis for polynomials of degree n is (f) k (t) =t k (0 < k < n). The 

Bernstein basis is <j> k (t) = Z Kn {t) = f ™ j (1 - t) n ~ k t k (0 < k < n). The 

Chebyshev basis is 4>k{t) = Tk(t) (0 < k < n), where Tk(t) is the Chebyshev 
polynomial of the first kind generated by the recurrence relation 





To(t) 
Ti(t) 
T k+1 (t) 



t 



2tT k (t) -r fc _i(t) for k > 1. 



1 



(9) 
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Another way to define the Chebyshev polynomials of the first kind is through 
the identity 

Tfc(cosa;) = cos/ca. (10) 

This second definition shows, in particular, that all zeros of Tk(t) lies in 
[—1,1]. It also shows that — 1 < Tk(t) < 1 for any —1 < t < 1. 

The rest of this article shows that the power, Bernstein, and Chebyshev 
bases all satisfy these basis properties. The values #'s of the three bases 
and their corresponding bounding polygons are also derived as these values 
dictate the efficiency of KTS operating on such bases. The upper bound of 
the efficiency of KTS is lowest when it operates on the basis with the smallest 

e. 



3.1 Bounding polygons 

The choices of I and h and the definitions of bounding polygons of the surface 
S = {/(«, v) : I < u, v < h}, where f(u, v) is represented by one of the three 
bases, that satisfy the required properties are as follows: For Bernstein basis, 
the convex hull of the coefficients (control points), call it P\, satisfies the 
requirements for / = and h — 1. The convex hull P\ can be described as 

Pi = < Ci i Si i '■ s « = - s « - 1 f • ( U ) 

I i,j i,j ) 

For power and Chebyshev bases, the bounding polygon 

Pz = S c 00 + ^2 CijSij : -1 < s i:j < 1 \ (12) 

I i+j>0 J 

satisfies the requirements for I — — 1 and h = 1. Note that P<i is a bounding 
polygon of S in the Chebyshev case since |T fc (t)| < 1 for any k > and any 
t E [—1, 1]. Determining whether G P 2 is done by solving a small linear 
programming problem. To determine if G Pi, the convex hull is constructed 
by conventional method and is tested if it contains the origin. 

The afline and translational invariance of P\ and P 2 for their respective 
bases can be verified as follows: Let 

m n 

g{u, v) = Af(u, v) + b = YY^ c'ijM^jiv)- 

i=0 j=0 
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For the Bernstein basis, by using the property that Y^k=o^k,n{t) — 1; it is 
seen that c^- = Ac^ + b for all c^-'s. Therefore, the bounding polygon of 
{g(u,v) : <u, v < 1} is 

p i' = |Z) t i a « : Z) a « = 1 ' ^ a «^ 1 } 

A ^2 °ij s ij + b Si i : = - - 1 f 

= {Ax + b : x e Pi} . 

For the power and the Chebyshev bases, note that 0o(w)0o(^) = 1 for 
both bases. Hence, c' 00 = Ac 00 + b and Cy = Acy for z + j > 0. The bounding 
polygon of {g{u, v) : < u, v < 1} for this case is 

^2 = <Uo + E c ^ s ^i : - 1 < < 1 f 
I *+i>o J 

{v4c o + b + ^ AcySij : -1 < < 1 > 
i+j>0 J 

= < A I c 00 + ^ cySij J + b : -1 < Sy < 1 
= {Ar + b: x E P 2 } ■ 

3.2 The size of the bounding polygons compared to 
the size of the bounded surface 

Item [2c] of the basis properties in effect ensures that the bounding polygons 
are not unboundedly larger than the actual surface itself lest the bounding 
polygons lose their usefulness. The value 6 also can be used as a measure 
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of the tightness of the bounding polygon. Recall from Theorem 12.11 that the 
efficiency of KTS depends on 9. 

Since the bounding polygons Pi and P 2 are defined by the coefficients 
of /, our approach to derive 9 is to first derive £, a function of m and n, 
satisfying 

||ctf|| < ^max_ ||/(«, v)\\ 



l<u,v<h 



for any coefficient Cjj of /. But the following lemma shows that one needs 
only derive the equivalent of £ for univariate polynomial to derive £ itself. 

Lemma 3.1. Assume there exists a function h(n) such that 

INI <%) max||<7(0|| (13) 

l<t<h 

for any hi (i = 0,1,..., n), and any univariate polynomial g(t) = 

XXtA'M*)- Then 

\\cij\\ < h(m)h(n) max ||/(m,i>)|| , (14) 

l<u,v<h 

for any Cij (i = 0,1, ... ,m; j = 0,1, ... ,n), and any bivariate polynomial 
f(u,v) = J2Zo J2]=o c ij4>i(u)(pj(v). 

Proof. Let f(u, v) = YliLo Ej=o c ij4 > i( u )4'j( v ) be an arbitrary bivariate poly- 
nomial. For any i = 0, 1, . . . , m, define g io (v) = YJj^^oAA 11 )- A PPly in g 
(USD to g io (v) yields 

llQojII < K n ) max \\g io (v)\\ , (15) 

l<v<n 

for any j = 0,1,..., n. Let v* = argmaxi<„< ft ||#i (f)||- Define l io (u) = 
E™oEJ=o c i3^(«)^K)- Applying (HH) to l io (v) yields 



j=0 



< h{m) max. \\l io (u)\\ , (16) 

Ku<h 



for any i = 0, 1, . . . , m. Consequently, by combining (fT5"|) and ffTBT) . 

||Qo,jII < K n ) maxK^ft \\g io (v)\\ = h(n) ||ft « )|| = h(n) E"=o ^o.A'K) ^ 
h(m)h(n)max.i< u < h \\l io (u)\\ < h(m)h(n)maxi< U;V < h \\f(u,v)\\. □ 
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We proceed to derive £ of the three bases, starting with the Bernstein 
basis. The following lemma regarding the product of two polynomials in 
Bernstein basis is needed to find £ for Bernstein case. 

Lemma 3.2. Let 

/(*) = Eto^nit), < t < 1 



and 
Then 

where 



n+n' 

f(t)g(t) = J2 b i Z i,n+n'(t), 
i=0 



\bi\ < max \ci\ ■ max \c[\ 

i i 

Proof. Straightforward arithmetic shows that 



n \ I n' 



min(n,i) ■ , ■ , . , 

/..,../ \ ' n-r)-i,- 

fc=max(0,i— n') 



n + n' 



Taking absolute value on both sides and bounding \ct\ (resp. |c^_ fc |) with 
maxj \ci\ (resp. max, \c' i \) gives 



min(n,i) 

< max |cj| • max |c-| 



n \ I n' 



fc=max(0,i— n') 



k J V i — k 

n + n' 

i 



Recall the combinatorial identity 



. min(ra,i) / \ / / 

n + n \ / n W n 

i = 2^ [ k [ i-k 

7 fc=max(0,i-7i') v 7 v 

Hence, the lemma follows. □ 
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With the above lemma, we are ready to derive £ of the Bernstein basis. 
Theorem 3.3. Let f(t) be a polynomial system 

/(*) = EtodZiAt), o < t < 1, 

where q G M d . The norm of the coefficients can be bounded by 



H| <£ B (n) t max ||/(f)||, (17) 



where 



Un) = J2 II = 

i=0 j=0,l,.,i-l,i+l,..,n J* 

Remark. An inequality in the other direction, namely, that 

max II f (t)\\ < max ||cj , 

is a well-known consequence of the convex hull property of Bernstein 
polynomials [2]. 

Proof. By definition of infinity norm, it suffices to prove the lemma for the 
case Cj G K. Therefore, it is assumed that d = 1 for the rest of this proof. 

Let tj = j/n (j = 0, 1, . . . ,n). Define a matrix A G R( n+1 ) x ( n+1 ) having 
element 

= Zi t n{tj)- 

Define the vectors c = (c , c 1 , . . . , c n ) T and / = (f(t ),f(t 1 ),...,f(t n )) T . 
Observe that 

Ac = f. (18) 

We claim that A is invertible. In particular, we show that the linear system 
Ax = b has solution for any arbitrary b G M™ +1 . Due to the definition of 
A, solving the system Ax = b is equivalent to finding the coefficients of the 
polynomial 

n 

git) = ^x t+l Z i>n (t) (19) 

t=0 
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with the property that g(t ) = 6 1( g{t\) — b 2 , ■ ■ ■ , g(t n ) = b n+ %. The polyno- 
mial g satisfying such property is the Lagrange interpolant 



'i+i 



j=0 



./' =a,...,j-i,j+i,...,n i 



n 



(20) 



Transforming (|20|) to the Bernstein basis yields the solution x. 

Knowing that A is invertible, we multiply both sides of (fl8|) by A^ 1 , 



(21) 



and hence, for any i = 0,1, 

\(k\ 



< \\c\\ 

< Wa-'I 



< lU^II ■ max 1/(01 • 

_ 11 11 (:0<K1 ' Wl 



(22) 



Comparing ( TTTI) to ( |22l) . it is seen that the final step is to show that \\A 1 1| < 
0>-i"). 

Observe that the ith column of A^ 1 is A 1 ej, where denotes the ith 
column of the identity matrix. Let gi(t) be a polynomial in the Bernstein 
basis and let {c' if } be its coefficients. With similar reasoning as the above, 





(c> \ 


= A- X 


/ ft(*o) ^ 


(23) 




\<) 




\ 9i(t n ) J 




But (123]) implies that the ith column of A 1 , A 1 e 


i, are the coefficients of g,i 


such that, for j — 0,1, . . 


,n, 










9i(tj 


Hi 


3 = h 
3 ^ i- 


(24) 



The following Lagrange interpolant gi satisfies ( |24l) : 
9i(t) = 



n 

j=0,...,i-l,i+l,-,n * ] 



n 

j'=0, l,i+l,...,n 



n ~ 3 1 3_ 



(25) 



11 



Note that each term of the product in (j2"51) is a polynomial in Bernstein basis 
with coefficients (n — j) j{% — j) and j / (i — j). Applying Lemma 13.21 to (1251) 
shows that 

M.i ii tt max{ I n — j I , I j I ) , , 

\\A~ 1 e i \\ < 11 ^ , h ulJ . (26) 

j. j. \i — j\ 

j=0,X,...,i—l,t+l,...,n 1 1 

Since (1261) holds for any column i of A -1 , the lemma follows. □ 

Next is the derivation of £ of the Chebyshev basis. The following identity 
is useful for this derivation: 

( i ± j 

^Ti^Tjfa) = I n i = j = (27) 
k=i (n/2 1=3^ 0, 

for i, j — 0, . . . ,n — 1, where (k = 1, 2, . . . , n) are the n zeros of T n (t). 
Theorem 3.4. Let f(t) be a polynomial system 

n 

f{t) = Y J c l rm, 

where q G M d . TTie norm of the coefficients can be bounded by 

IM|<V5 t m^ i ||/(t)||. (28) 

Proof. By definition of infinity norm, it suffices to prove the lemma for the 
case q G R. Therefore, it is assumed that d = 1 for the rest of this proof. 

Let tj (j = 1, 2, . . . , n + 1) be the n + 1 zeros of T n+ i(t), which lie in 
[—1,1]. Define a matrix A G R( n+1 ) x ( n+1 ) having element 

Define the vectors c = (c , c x , . . . , c n ) T and / = (f(t ),f(t 1 ),...J(t n )) T . 
Observe that 

Ac = f. (29) 

By fl27j), 

v4 T v4 = diag (n + 1, (n + l)/2, (n + l)/2, . . . , (n + l)/2) , 
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which implies that A is invertible and 

A~ l A~ T = diag (l/(n + 1), 2/(n + 1), 2/(n + 1), . . . , 2/(n + 1)) . (30) 
The equation (1311 implies 



Lr 1 !! . = VV(n + i)- (si) 



Finally, from (|29l) and (|3T|), 

| Q | ^ 1 1 ^ 1 1 2 

< IU" 1 



2 

< Vn + 1 U _1 | 



< y/n + 1 |U _1 |L max 

11 112 t:-KKl 



max \f(t)\. □ 



fc-KKl 



Last is the power basis. Our approach to derive £ of power basis is to 
derive the relationship between the coefficients of a polynomial in power basis 
and the coefficients of the same polynomial but written in Chebyshev basis. 
By using this relationship and Theorem 13.41 £ of the power basis can be 
computed. 

Lemma 3.5. Let f be a univariate polynomial such that 

n n 
i=Q i=0 

In other words, {aj} are the coefficients of f when written in the power basis 
and {cj} are the coefficients of f when written in the Chebyshev basis. Then 

3 n+1 - 1 
| ai| < max \cj\, 

2 j=0,...,n 

for any i — 0, . . . , n. 

Proof. Let D = [dij] be the n + 1-by-n + 1 matrix such that 

a = Dc, 
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where a = (a , a%, . . . , a n ) T and c = (c , Ci, . . . , c n ) T . Note that 

j 

^(*) = y^j+ij+i^- 

i=0 

Recall the recurrence relation Tj(t) = 2t2}_i(t) — Tj_ 2 (t). It follows from 
this recurrence that 

K+ij+il < 3 J '. (32) 

That is, when Tj(t) is written in the power basis, the resulting coefficients (of 
power basis) is less than or equal to 3 J . The inequality (I3"2"j) can be verified 
by induction on the recurrence relation. Since the entries in the (j + l)th 
column of D is bounded by 3 J , we have, from geometric sum, 

\\D\\ < (3 n+1 - l)/2. 

The lemma follows from ||a|| < \\D\\ \\c\\. □ 

Theorem 3.6. Let f(t) be a polynomial system 

n 
i=0 

where q G The norm of the coefficients can be bounded by 

on+l i 

w *-^-.3&. |l/wl1 - (33) 

Proof. Follow directly from Theorem 13.41 and Lemma 13.51 □ 

Having ^ for each of the three bases, the values of 9 for the three bases 
can now be derived. 

Corollary 3.7. Let 

m n 

f(u,v) = ^^CijZi^ujZj^v), 

i=0 j=0 

where Cy G IR 2 (i = 0, 1, . . . , m; j = 0, 1, . . . , n) . Let Pi be the convex hull of 
{cij}. Then, for any y G Pi, 

IMI < £s(W)G?(n) max \\f(u,v)\\ . 

0<u,v<l 
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Proof. By the convex hull property of Bernstein polynomials, \\y\\ < 
maxjj || Cij || for any y G Pi. The corollary then follows from Theorem 13.31 
and Lemma 13. 11 □ 



Corollary 3.8. Let 



/(«,«) = ££ 



i=0 j=0 

where G M 2 (z = 0, 1, . . . , m;j = 0, 1, . . . , n). Let 



I c 00 + ' 

I i+j>0 



Then, for any y G Pi, 



I || ^ (m+l)(n + l)(3 m+1 -l)(3» +1 -l) 

Z — l<u,v<x 



m n 



Proof. For any ?/ G P 2 , 

i=0 j=0 

The corollary then follows from Theorem 13.61 and Lemma 13. 11 □ 
Corollary 3.9. Let 

m n 
i=0 j=0 

where Gt 2 (i — 0, 1, . . . , m; j — 0, 1, . . . , n). Let 
Pi 



\Cqo+ CySij : -1 < Sij < 1 > . 
I j+i>o J 



T/ien, /or an?/ y E P2, 



\\y\\ < 2(m + l)(n + 1) max u)|| . 

— 1<M,^<1 

Proof. For any w G P2, 

m n 
i=0 j=0 

The corollary then follows from Theorem 13.41 and Lemma 13. 11 □ 
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4 Relationship between the bounding poly- 
gon of the power basis and that of Cheby- 
shev basis 

Let P$ denote the bounding polygon P 2 computed from the power basis rep- 
resentation of a polynomial and P% denote Pi computed from the Chebyshev 
basis representation of it. The results from previous section show that the 
value 9 of P 2 C is smaller than 9 of Pf. This only implies that the worst case 
of P| is better than the worst case of Pf. Comparing the values of #'s of 
the two does not indicate that P 2 C is always a better choice than Pf for every 
polynomial. The following results show, however, that P 2 C is, in fact, always 
a better choice than P 2 . Specifically, this section shows that for any given 
polynomial, its bounding polygon P 2 C is a subset of its bounding polygon P|\ 
The following two lemmas show that when representing monomials t k in 
Chebyshev basis, each coefficient is nonnegative, and the sum of all coeffi- 
cients are exactly 1. These results are useful in relating Pf to P 2 C . 

Lemma 4.1. Let d^ 's (k = 0, 1, . . .; i = 0, 1, . . . , k) be the numbers satisfying 
t k = Ei= dkiTi(t). Then 

dki > 0, 

for any k = 0, 1, . . . and any i = 0, 1, . . . , k. 

Proof. We prove the lemma by induction on k. The base cases k = and 
k = 1 are trivial. For the inductive step, for any k > 1, 

t k+1 = t-t k 

k 

= t^2d ki Ti(t) 

= ( 2tT ^) - r «-i(0) + E Y Tl - l{t) + 4otT °W- 

i=l i=l 

By © and noting that tT (t) =t = T 1 {t), 

t k+1 = ty^W + ty^W + «W 

8=1 8=1 

= E nr r '(*) + E t dk,i+l Ti{t) + (t + 4o ) Tl{t) + T Tl 

i=k 8=2 ^ ' 
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Hence, 

4,i-i/2, i = k,k + l, 

, _ j (dk,i-i + dk,i+i) /2, i = 2, . . . , k — 1, , , 

dfc+w_ i 4 2 /2 + 4o, i = l, (34) 

4i/2, i = 0. 

But since d&, > for any i — 0, . . . , k by the induction hypothesis, f[54"|) shows 
that > for any % = 0, . . . , k + 1. □ 

Lemma 4.2. Lei 's ('A; = 0, 1, . . .; i = 0, 1, . . . , k) be the numbers satisfying 
t k = J2i = o d kiTi{t). Then 

k 

dki = i, 

i=0 

/or any fc = 0, 1, . . . and any i — 0, 1, . . . , k. 

Proof. We prove the lemma by induction on k. The base cases k = and 
k = 1 are trivial. For the inductive step, the same reasoning as in the proof 
of Lemma [4.11 shows that dk+ij is as in (|34"|) for any k > 1. Therefore, 

E v-^ dfc.i-i \ d& j_i + dfe i+i / dfc2 , , \ , dki 
4 +M = + L - 2 + (l- + dfe0 ) + T 

i=0 i=fc «=2 v 7 

fc 



dfcj — 1, 



i=0 

by the induction hypothesis. □ 

Theorem 4.3. Let / : R 2 ^ R 2 k a bivariate polynomial. Its bounding 
polygon P| is a subset of its bounding polygon Pf . 

Proo/. Let / = /(u,v) = ELoELo^^' = J2T=oT^=o c ij T i( u ) T j( v )^ 
where a^i € M n (A; = 0, 1, . . . , m; / = 0, 1, . . . , n) are the coefficients of / when 
written in the power basis and G W 1 (i = 0,1, ... ,m; j = 0,1, ... ,n) 
are the coefficients of / when written in the Chebyshev basis. Let d^'s 
(k = 0, 1, . . .; i = 0,1, ... ,k) be the numbers satisfying t k = Yli=o dkiTi(t). 
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Hence, 



f(u,v) = akl [ d ^ T i( u ) J I d ii T i( v ) 

k=0 1=0 \i=0 / \j=0 



m n m n 



i=0 j=0 k=i l=j 

Therefore, = J^fcLi S™=j a kidkidij. This means that P 2 C can be written as 

{m n m n 

a 00 d 00 d 00 + a-kidkodio + dkidkidijSij : — 1 < < 1 ^ , 

fc=l Z=l i+j>0 fc=i «=j 

But since rf o — 1) 



m n n m n 



P 2 C — ^ a o + dkidkodio + o J kidkodijS j+ 

k=l 1=1 j=l k=0 l=j 

m m n m n m n 

a kidkidioS i0 + 2j a kidkidijSij : -1 < s^- < 1 

i=l fc=i i=0 i=l j=l fc=i i=j 

m n m n I 

a °° + akldk ° dl ° + a kl d kodljSoj + 

k=l 1=1 k=0 1=1 j=l 

m n k m n k I 

YjYjYj a kidkidioS i0 + "kidkid/j*;, : -1 < < 1 

k=l 1=0 i=l k=l 1=1 i=l j=l 

n I m k 

a °° a ° l dl 3 S °j + afc0 5-/ dkiSi0 ^~ 

1=1 j=l k=l i=l 

m n / I k k I 

a fc« I dkodi + d k0 (i/jSoj + c?io G^io + cfe 

fe=l i=l \ j'=l i=l i=l j'=l 

-1 < sn < l\. 

By Lemma H~T1 and Lemma H~2| it is seen that —1 < Yuj=i dij$oj, Yli=i dki$io < 
1. In addition, using the fact that < 1, for any i = 0, . . . , m and any 
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j = 0, . . . ,n, together with Lemma [4.11 and Lemma 14. 2\ it is seen that 
dkodio + d k0 Y^j=i dijSoj + <^o Yli=i dkiSio + Yli=i dki Y^j=i < |c4o^zo| + 

\dw\ Y^j=i \dij\ + \dio\ Yli=i \dki\+Yli=i \dki\ X^j=i \dij\ — 4o^o+4o Ej=i^i+ 
dio J2i=i dki + J2i=i dki Y!j=i dij = \ J2i=o d ki) (Ei=o^i) = L Therefore, 

4.1 Reparametrization 

The last nontrivial basis property that warrants detailed discussion is the 
issue of efficient reparametrization. Reparametrizing polynomials in power 
basis is straightforward from the binomial theorem. Polynomials in other 
bases, on the other hand, may not be as simple to reparametrize. The details 
of the process for polynomials in Bernstein and Chebyshev bases are covered 
in this section. 

4.1.1 Reparametrization of polynomials in Bernstein basis 

There is more than one algorithm to compute the reparametrization with 
[0, l] 2 of a bivariate polynomial in Bernstein basis. We describe one method 
here. Our method makes use of a program that, given a^s, c, d, e, /, g, h, 
k, and I, computes /3y's satisfying 



m n 



E E a ^°y + + /r~V + + ir~ j = E E W**- 

Such conversion can be done in O ((mn) 2 ) by generalizing Horner's rule. We 
leave the details of the conversion to the reader. Let X denote {(u,v) : 
u° — r < u < u° + r, v — r < v < v + r}. To compute the coefficients {%} 
°f {fx{u, v) : < u, v < 1}, the [0, l] 2 -reparametrized surface of {f(u,v) : 
u° — r < u < u° + r, v — r < v < v + r}, first substitute u = 2ru + u° — r 
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and v = 2rv + v — r into /, yielding 

i=0 j=0 ^ / \ J J 



= E E f T ) f • ) c -( 2 ™ + m ° - ^ C 1 - ( 2 ™ + m ° - ^)) m "* ■ 

i=o j=o ^ / \ J / 

(2rv + v° - r) j (l - (2rv + v° - r))"^ 

= EE(?) (")c,,((. + r). + (n -r)(l-n)) J . 
i=o j=0 ^ J \ J J 

((1 - u° - r)u + (1 - u° + r)(l - M)) m "' • 

((v° + r)v + (v° - r)(l - v)) j ■ 

((1 -v°- r)v + (1 - v° + r)(l - . (35) 

Substituting u — u/{u + 1) and t) = v/(v + 1) into (155]) yields 

= t ^YY( m V n )c ij ((u + r)u + u°-r) i - 



((l-u -r)n + l-M° + r) m_i - 
((u° + r)u + w° - r) J • 
((1 - v° - r)v + 1 - v° + r)" _i . (36) 

n 

EE^^> ( 3? ) 



(Hl) m (Hl)^ n n 

v ' v 7 2=0 J=0 

where (1371) is obtained from (1361) by the conversion program mentioned above. 
Substituting u and v back into (1371) to see that 

m n 

/(«,«) = E !>««'(! ~ ^) m ^'( 1 - (38) 
i=o j=0 

Therefore, Cy = r y i j/(C(m,i)C(n,j)) are the control points of /x where 
C(m,i) 
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4.1.2 Reparametrization of polynomials in Chebyshev basis 

Let a, b, d and e be scalar constants. The reparametrization with [—1, l] 2 of 
a bivariate polynomial in Chebyshev basis can be computed if the values of 
Aifc's (i — 0, 1, . . . , m) satisfying 

i 

T i {at + b) = Y,\kT k {t) 

k=0 

and the values of fij k s (j = 0, 1, . . . , n) satisfying 

j 

Tj(dt + e) = ^2/j jk T k (t) 

k=0 

are known. Note that 

i j 

Ti(au + b)Tj(dv + e) = ^ ^ \ ik ^ jkl T k (u)T k ,(v), 

k =0 k'=0 

which is adequate to find the reparametrization. The values of a, b, d, and e 
are determined by the iif-domain of the surface to be reparametrized. 
To compute Xi k s, observe that for i > 1, by fl9j), 

T i+1 {at + b) = 2(at + b)Ti{at + b) -T^at + b) 

i i — 1 



2(at + b)J2 ^ikT k {t) - Vi,fc7fc(t) 

fc=0 fc=0 
i i i— 1 

^ 2aA ife tT fc (t) + £ 26A ifc T fc (t) - ^ A,_ life T fc (t) 

k=0 k=0 k=0 

% 

2a\ i0 tT {t) + aX ik (2<T fc (t) - T fc _i(t)) + 
fc=i 

i— 1 i i—1 

^2 a\ k +iT k {t) + ^ 2b\ ik T k (t) - ^ \-i,kT k {t) 

k=0 k=0 k=0 

i 

2aX i0 T 1 (t) + y^ j aX ik T k+1 (t) + 

k=l 

i—1 i i—1 

Y, aX ^+iT k (t) + Y,2b\ ik T k {t) - \-iMt), (39) 



k=0 k=0 k=0 
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and 



T {at + b) 
T^at + b) 



T (t), 

aT\(t) + 6T (t). 



(40) 
(41) 



The equalities (139|) . (140|) . and (j4"Tj) yield a recurrence relation of Ajfc's that 
can be used to compute their values. The values of fijk's can be computed 
similarly. 

5 Computational results 

Three versions of KTS algorithms are implemented in Matlab; one operat- 
ing on the polynomials in power basis, one on Bernstein basis, and one on 
Chebyshev basis. They are tested against a number of problem instances 
with varying condition numbers. Most of the test problems are created by 
using the normally distributed random numbers as the coefficients Cjj's of 
/ in Chebyshev basis. For some of the test problems especially those with 
high condition number, some coefficients are manually entered. The resulting 
Chebyshev polynomial system is then transformed to the equivalent system 
in the power and the Bernstein bases. Hence the three versions of KTS solve 
the same polynomial system and the efficiency of the three are compared. 
The degrees of the test polynomials are between biquadratic and biquartic. 

For the experiment, we use the algorithm by Jonsson and Vavasis [3] to 
compute the complex zeros required to estimate the condition number. Table 
CD compares the efficiency of the three versions of KTS for the test problems 
with differing condition numbers. The total number of subpatches examined 
by KTS during the entire computation and the width of the smallest patch 
among those examined are reported. The results do not show any one version 
to be particularly more efficient than the others although the Chebyshev basis 
has better theoretical bound than the other two. 

Since the types of test polynomials may affect the relative efficiency of the 
three versions of KTS, another experiment is performed on degree 6 univari- 
ate polynomials generated by different methods. Since Section [5] shows that 
the Chebyshev basis always gives tighter bounding polygons than the power 
basis, this experiment only compares between the Chebyshev and Bernstein 
bases. Table [2] and Table [3] show the results of this experiment. The poly- 
nomials are generated as follows. The "rand" polynomials are generated by 
interpolating points whose x-coordinates are evenly spaced between —1 and 
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Table 1: Comparison of the efficiency of KTS algorithm operating on the 
power, the Bernstein, and the Chebyshev bases. The number of patches 
examined during the course of the algorithm and the width of the smallest 
patch examined are shown for each version of KTS. 





Power basis 


Bernstein basis 


Chebyshev basis 


cond(/) 


Num. of 


Smallest 


Num. of 


Smallest 


Num. of 


Smallest 




patches 


width 


patches 


width 


patches 


width 


3.8 x 10 3 


29 


.125 


17 


.0625 


21 


.125 


1.3 x 10 4 


13 


.125 


17 


.0625 


13 


.125 


2.5 x 10 5 


49 


.0625 


21 


.0625 


45 


.0625 


1.1 x 10 6 


97 


.0313 


65 


.0313 


85 


.0313 


3.9 x 10 7 


89 


.0313 


81 


.0313 


89 


.0313 



1, inclusive, and whose y-coordinates are normally distributed random num- 
bers. The "sin" ones are interpolations of sm(ax + b) at evenly spaced points 
between —1 and 1, inclusive, where a and b are normally distributed random 
numbers. The "sin-L" ones are the same as the "sin" ones except that least- 
squares interpolation is used instead. The "sinw" (resp. "sinw-L") ones are 
generated in the same way as the "sin" (resp. "sin-L") ones but with the 
function sin(6ax + b). Table [2] compares the number of test polynomials of 
each type where one basis yields tighter bounding intervals than the other. 
Table [3] shows the number of test polynomials of each type that bounding 
intervals of each basis have at least one endpoint exactly at the boundary of 
the ranges of the polynomials. The results show that the Chebyshev basis 
is decidedly better for "rand" polynomials, is about the same for "sin" ones, 
but is worse for the rests of the polynomials than the Bernstein basis. 

6 Conclusion 

Three common bases, the power, the Bernstein, and the Chebyshev bases, 
are shown to satisfy the required properties for KTS to perform efficiently 
In particular, the values of 9 for the three bases are derived. These values 
are used to calculate the time complexity of KTS when that basis is used to 
represent the polynomial system. The Chebyshev basis has the smallest 9 
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Table 2: The numbers of test polynomials out of 1000 that bounding intervals 
associated with the Bernstein basis is tighter than the those associated with 
the Chebyshev basis, and vice versa. 



Poly, type 


Num. that Bernstein is tighter 


Num. that Chebyshev is tighter 


rand 


1 


999 


sin 


963 


37 


sin-L 


960 


40 


sinw 


436 


564 


sinw-L 


998 


2 



Table 3: The numbers of test polynomials out of 1000 that bounding intervals 
associated with the Bernstein basis and those associated with the Chebyshev 
basis having at least one endpoint exactly at the boundary of the ranges of 
the polynomials. 



Poly, type 


Num. Bernstein with 
exact endpoint 


Num. Chebyshev with 
exact endpoint 


rand 


2 


13 


sin 


965 





sin-L 


972 





sinw 


330 





sinw-L 


999 


658 
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among the three, which shows that using KTS with the Chebyshev basis has 
the smallest worst-case time complexity. The computational results, how- 
ever, show no significant differences between the performances of the three 
versions of KTS operating on the three bases. It appears that, in average 
case, choosing any of the three bases do not greatly affect the efficiency of 
KTS. The experiment on univariate polynomials show that the Bernstein 
basis is more suitable for certain types of polynomials while the Chebyshev 
basis is better suited for other types. 
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